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Abstract 

In this work, the modeling of triaxially braided composites was explored through a semi-analytical 
discretization. Lour unique subcells, each approximated by a “mosaic” stacking of unidirectional 
composite plies, were modeled through the use of layered-shell elements within the explicit finite element 
code LS-DYNA. Two subcell discretizations were investigated: a model explicitly capturing pure matrix 
regions, and a novel model which absorbed pure matrix pockets into neighboring tow plies. The in-plane 
stiffness properties of both models, computed using bottom-up micromechanics, correlated well to 
experimental data. The absorbed matrix model, however, was found to best capture out-of-plane flexural 
properties by comparing numerical simulations of the out-of-plane displacements from single-ply tension 
tests to experimental full field data. This strong correlation of out-of-plane characteristics supports the 
current modeling approach as a viable candidate for future work involving impact simulations. 

1.0 Introduction 

Triaxially braided composites have been implemented into fan containment systems for jet engines due 
to the high-specific stiffness and strength of fiber-reinforced polymer matrix systems and the improved 
impact damage tolerance and delamination resistance of braided architectures. These types of composites, 
however, present significant disadvantages for numerical modeling and simulation. Lirst, previous impact 
studies have shown architectural dependence on damage evolution and damage patterns during a ballistic 
impact event (Ref. 1). Consequently, an appropriate modeling strategy must account for the non- 
homogeneity present at the mesoscale: the scale of the braided pattern. Second, standardized coupon testing 
used for metals and simpler composite architectures have been shown to be relatively ineffective for 
determining braided composite strength properties (Ref. 2). Thus, accurate strength properties required for 
standard finite element composite damage models are currently unavailable. Additionally, the modeling 
scheme must be computationally efficient in order to be suitable for modeling of dynamic events such as a 
fan-blade-out event. The current work presented proposes to improve on already developed semi-analytical 
methods utilizing layered-shell elements and existing composite material models within the explicit finite 
element code LS-DYNA (Ref. 3). The overall goal is to develop a modeling algorithm which can best 
capture mechanisms of deformation and damage unique to the braided reinforcement architecture. 

Analytical methods for determining the elastic properties of braided composites utilizing 
micromechanics were first developed as an extension of previous literature on woven composites. These 
typically involved first identifying the repeating unit cell, RUC, of the braided mesostructure. Ishikawa 
and Chou developed mosaic, fiber undulation and bridging models, modifying classical laminate theory to 
account for fiber undulations (Refs. 4 and 5). Naik and Shembekar (Refs. 6 and 7) extended the one- 
dimensional mosaic model into two dimensions, whereas Yang et al. (Ref. 8) improved on accounting for 
fiber undulation in a fiber inclination model. Branch et al. (Ref. 9) extended a fiber inclination model into 
three dimensions. Redman and Douglas (Ref. 1 0) implemented the mosaic model using the rule of 
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mixtures and classical lamination plate theory (CLPT). More advanced discretizations of the RUC were 
developed by Ivanov and Tabei (Refs. 1 1 to 13), splitting it into four unique subcells. Huang (Refs. 14 to 
16) developed a micromechanical model modeling individual subcells as regions of two yams and two 
pure matrix regions and treating each as a unidirectional (UD) composite. Quek et al. (Ref. 17) developed 
an analytical model which consisted of one axial tow, two braider tows, and one matrix layer. Concentric 
cylinder modeling was used to determine the homogeneous properties of the RUC. Carey et al. (Ref. 18) 
divided the braided RUC into 13 distinct regions and utilized a volume weighted stiffness matrix to 
homogenize the properties of the RUC. It was noted in Masters et al. (Ref. 19) that simple analytical 
laminate models performed well for stiffness predictions, whereas Falzon (Ref. 20) noted that these 
representations were unable to accurately capture out-of-plane properties. 

There are three classifications of numerical methods for predicting the response of triaxially braided 
composites. The first includes those which model the RUC explicitly through finite element (FE) models 
(Refs. 21 to 23). In these approaches, the explicit models are used to determine the effective properties 
(stiffness/strength) of the RUC for use in structural simulations. A review of these types of numerical 
methods can be found in (Ref. 24). A second branch utilizes multiscale methodologies to li nk the explicit 
FE representation of the RUC with homogenized elements on the structural level (Ref. 25). Lastly, a 
semi-analytical method was developed (Refs. 26 and 27) based on a braiding-through-the-thickness 
method introduced by Cheng and Binienda (Ref. 28) utilizing layered shell elements, subcell 
discretizations, and the inherent CLPT implementation available in FE packages. Xiao et al. (Ref. 29) 
utilized a similar layered shell methodology with a unique subcell discretization. 

In the current work, the semi-analytical methodology is utilized and builds upon limitations in the 
previous models by Cheng and Binienda (Ref. 28) and Xiao et al. (Ref. 29). Section 2.0 will discuss the 
general methodology in the semi-analytical approach using subcell models and introduce the novel 
absorbed matrix model. Section 3.0 will cover the numerical implementation of the modeling procedure 
and Section 4.0 will verify the in-plane elastic properties of the subcell models. Finally, the efficacy of the 
absorbed matrix model will be demonstrated with novel simulations of single -ply transverse tension tests 
in Section 5.0. 


2.0 Subcell Modeling 

2.1 Previous Works 

The material system of interest is a [0°/+60°/-60°] triaxially braided composite made of 24 k carbon 
fiber bundles in the 0° tow direction and 12 k carbon fiber bundles in the ±60° braider directions. Braider 
geometries were formed by A&P Technologies and RTM molding was done by North Coast Composites 
(Refs. 1, 2, and 30). Fiber and matrix constituent sources and properties are provided in Section 3.1. A 
typical RUC, the smallest characteristic volume capturing the effective properties of this braided 
composite, is shown in Figure 1(a). The unit cell can be decomposed into 4 subcells as outlined in 
Figure 1(b). This partitioning corresponds to the existence of regions where axial (0°) and braider 
(+60°/-60°) tows are both present and regions where only the braider tows are present. Figure 1(b) only 
displays the fiber tows in the braided composite and not the pure matrix regions which exist between and 
around fiber tows. The boundary of subcells A and C are drawn to lie at the edges of axial fiber tows. The 
individual subcells can be subsequently discretized in a “mosaic” model of UD composite plies to capture 
the orientation and response of the fiber tows within the triaxial braid. Previous subcell modeling 
(Refs. 26 to 28, and 30) developed in collaboration with the NASA Glenn Research Center (GRC) and the 
University of Akron, shown in Figure 2(a), utilized subcells containing two unique UD plies: one 
containing the axial tow, and another containing the braider tow and surrounding matrix material. The 
two UD plies were assigned unique fiber volume fractions. These volume fractions were determined by 
matching global fiber volume fractions and assuming 80 percent fiber volume fraction in the 0° plies. The 
fiber volume fraction of the axial tows in this braided composite have been previously assumed to be 
80 percent (Ref. 26). 
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(b) ■0° H+60° B-60° 

Figure 1. — (a) Unit cell of the braided composite, (b) The subcell 
segments with braid angles shown. 




Figure 2. — (a) Subcell UD discretization of previous layered shell 
models developed at GRC and U of Akron (Refs. 26 to 28, 30). 
(b) The subcell model developed by Xiao et al. (Ref. 29). 


There were two limitations to these previous models which are addressed in this work. First, the latest 
revision of the semi-analytical model homogenized each subcell through the thickness, resulting in a 
single layer for each subcell (Ref. 27). Although the model was still able to capture in-plane stiffness and 
strength properties of the braided composite, unique flexural properties of individual subcells were lost, 
such as the unique tension-twist coupling of the unsymmetric subcell laminate documented by Kohlman 
(Ref. 2). Secondly, the fiber volume fraction of the ±60° UD plies in the model in Figure 2(a) were 
determined by meeting the global composite fiber volume fraction, but did not reflect the actual fiber 
volume fraction contained in the given subcells. 

The subcell model developed by Xiao et al. (Ref. 29) overcame these limitations by retaining the 
layered shell formulation for each individual subcell as well as developing a simple straight line model for 
determining ply thicknesses and UD fiber volume fractions. Xiao et al. also modeled fiber tows as unique 
UD plies, where the fiber volume fraction of each ply was identical to that of the fiber tow (~80 percent) 
and pure matrix pockets were modeled explicitly using the resin properties. Figure 2(b) displays the 
discretization, which hereon is referred to as the pure matrix model. Another significant advantage to the 
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Figure 4. — Cross-section micrograph of the braided 

Figure 3— Schematic highlighting the equivalence of composite, at lOOx magnification, displaying matrix 

+/-60 0 and -/+60° tows in subcells B and D. pockets adjacent and parallel to the fiber tows. 

discretization of the pure matrix model was that the B and D subcells containing only the braider tow 
plies were modeled as symmetric [-60°/+60°/-60°] and [+60°/-60°/+60°] laminates. Figure 3 shows a 
detailed schematic highlighting the overlapping of +60° and -60° tows in subcells B and D. Although the 
orientation of the overlap differs between the subcells, they both contain equal amounts of ±60° and 
-/+ 60° tow regions within the entire subcell. Although the local green/yellow triangular regions on the 
braid are locally unsymmetric, one can approximate that these two subcells are symmetric with respect to 
the subcell discretization. 

A limitation of the pure matrix model, however, is the fact that the pure matrix pockets contained 
between the UD plies do not accurately represent the true interaction of pure resin pockets with adjacent 
fiber tows. Figure 4 presents a cross-sectional micrograph of the braided composite with a red arrow 
indicating a pure matrix region adjacent to an axial fiber tow and a green arrow indicating matrix rich 
regions between axial and braider tows. From the micrograph images, more pure matrix regions are found 
lying in-plane to either the braider or axial tow rather than sandwiched between braider and axial tows as 
is modeled in the pure matrix model in Figure 2(b). Despite this deficiency, the pure matrix model is used 
as a comparative reference to a new absorbed matrix model. 

2.2 Absorbed Matrix Model 

In an attempt to combine aspects of both the previous semi-analytical models with that of the pure 
matrix model, it was decided that a suitable subcell discretization would contain the following elements: 

1 . Explicit modeling of UD layers 

2. Accurate subcell fiber volume fractions 

3. Symmetric B/D subcells 

4. Better matrix representation 

The absorbed matrix model utilizes the straight line model used by Xiao et al. for determining fiber 
volumes in each subcell for the axial and braider directions. Rather than modeling pure matrix pockets as 
unique layers which are in parallel (that is, with respect to in-plane loading) to the UD fiber tow plies, the 
volume of pure matrix is “absorbed” by the surrounding UD fiber tow plies to create a homogenized 
representation of the tow/matrix interaction. This representation of the tow/matrix interaction is notably 
similar to the original subcell model in Figure 2(a), however, the fiber volume fraction of the braider plies 
in this model are now thought to best capture the amount of fibers in the respective subcells. Figure 5 
compares the fiber volume fractions between the previous semi-analytical subcell models and those of the 
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pure matrix and absorbed matrix models. Both the pure and absorbed matrix models have similar volume 
fractions since they are derived based on the same straight line model approach. Figure 5 shows fiber 
volume fraction totals (blue) for subcells A and B and fiber volume fractions segregated by direction (red 
and green). It should be noted that the volume fractions shown for the axial and braider plies are given 
with respect to the volume of the subcell, thus are not ply-level fiber volume fractions discussed later. 

The geometric parameters required to use the straight line model for this braided composite system 
are given in Table 1. Both the pure matrix model and the absorbed matrix model are applied in the 
sections to follow for testing both in-plane stiffness properties and single ply coupon tests. The resulting 
subcell discretizations are presented in Tables 2 and 3, respectively. The necessary equations and figures 
for the straight line model used for computing the thicknesses in Tables 2 and 3 are presented in the 
Appendix. The thicknesses in Tables 2 and 3 are with respect to a total shell thickness of unity. It should 
be noted that the absorbed matrix model assumes that pure matrix pockets in subcells A and C are 
assimilated only into the ±60° plies. Thus, the axial UD ply in subcells A and C has a fiber volume 
fraction of 80 percent for both the pure matrix and absorbed matrix models. The next section will cover 
the numerical implementation of both models within LS-DYNA. 


Fiber volumes for subcells A and B 
H Subcell-overall 



Previous subcell models (Ref. 27) Pure matrix/absorbed matrix model 

Figure 5. — Comparison of fiber volume fractions at the subcell level between 
the previous subcell models and the pure matrix and absorbed matrix models. 


TABLE 1.— GEOMETRICAL PARAMETERS USED FOR 
CALCULATING THE LAMINA THICKNESSES 
(T700/PR520 SYSTEM) 


Label 

Description 

Value 

W 

Width of RUC (m) 

0.0089 

W A 

Width of cell A (m) 

0.004201 

w B 

Width of cell B (m) 

0.004765 

h 

Ply thickness (m) 

0.00056 

Vf.UD 

Fiber volume (tow) 

0.8 

n a 

Number of fibers in axial tow 

24000 

n b 

Number of fibers in braider tow 

12000 

d a 

Diameter of fiber filament in axial tow (m) 

0.000007 

db 

Diameter of fiber filament in the braider tow (m) 

0.000007 

L 

Height of unit cell (length) (m) 

0.0051 

0 

Braid angle (degrees) 

60 
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TABLE 2.— SUBCELL A & B DISCRETIZATION FOR 
THE PURE MATRIX MODEL 


Subcell A Lay-up 

Angle, 

degree 

Vfi 

percent 

Thickness 

Braider tow 

-60 

80 

0.2334 

Matrix 

— 


0.02123 

Axial tow 

0 

80 

0.49074 

Matrix 

— 

— 

0.02123 

Braider tow 

60 

80 

0.2334 

Subcell B Lay-up 

— 

— 


Braider tow 

-60 

80 

0.1174 

Matrix 

— 

— 

0.2652 

Braider tow 

60 

80 

0.2348 

Matrix 

— 

— 

0.2652 

Braider tow 

-60 

80 

0.1174 


TABLE 3.— SUBCELL A AND B DISCRETIZATION FOR 


THE ABSORBED MATRIX MODEL 


Subcell A Lay-up 

Angle, 

degrees 

Vf 

percent 

Thickness 

Braider tow 

-60 

73.0 

0.255 

Axial tow 

0 

80.0 

0.49 

Braider tow 

60 

73.0 

0.255 

Subcell B Lay-up 

— 



Braider tow 

-60 

37.5 

0.25 

Braider tow 

60 

37.5 

0.50 

Braider tow 

-60 

37.5 

0.25 


3.0 Numerical Implementation 

3.1 Material Models 

The work presented here focuses on single-ply coupon tests which remained within the elastic range 
of the composite, thus all plies in both subcell models required only elastic properties. However, in 
preparation for future work to identify a top-down approach for determining ply-level strength, viable 
candidates for material models were investigated and utilized in the work. The matrix material was 
assumed to be elastic, perfectly-plastic and was implemented within LS-DYNA using material type 24, 
*MAT_PIECEWISE_LINEAR_PLASTICITY. The matrix properties were taken from Table 4. The 
efficacy of the elastic, perfectly-plastic assumption was not validated here, and it should be noted that the 
explicit modeling of the matrix was not necessary for the absorbed matrix model. Two resident material 
models of LS-DYNA were identified as candidates for modeling the unidirectional composite material: 
Mat_54, a linear-elastic composite damage model, and Mat_58, a non-linear composite model utilizing a 
continuum damage mechanics (CDM) approach. The former utilizes the Chang-Chang failure criterion 
and models elastic/brittle failure in the fiber direction and an elastic, perfectly-plastic response in the 
transverse and shear matrix directions. The latter model uses a modified-Hashin failure criterion and a 
strain-driven CDM formulation developed by Matzenmiller (Ref. 31). A purely elastic composite 
response was simulated by setting arbitrarily high strength values in Mat_58 for all directions in the 
composite damage model, thereby requiring only updated elastic properties. For these stiffness values, a 
bottom-up micromechanics approach was employed to calculate the moduli and Poisson’s ratio of the 
unidirectional layers in each of the subcells. In this way, the effective properties of the UD ply were 
determined based on the constituent properties of the fiber and matrix. The micromechanics software 
MAC/GMC developed at NASA GRC (Ref. 32) was used to compute the effective properties (elastic 
moduli, Poisson ratio) of the UD lamina based on the constituent properties provided in Table 4. It was 
hypothesized to be sufficient to use this bottom-up approach since the UD discretization of the sub-cell 
accurately captured the fiber contributions in all material axes present in the braided composite 
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TABLE 4.— CONSTITUENT MATERIAL PROPERTIES 
(REFS. 26 TO 28 AND 33) 


Material 

Density, 

g/cm 3 

Ell, 

GPa 

E22, 

GPa 

Poisson, 

vl2 

G12, 

GPa 

T700 

1.8 

230 

15 

0.2 

27 

PR520 

1.256 

4 

4 

0.4 

1.47 

5208 

1.263 

3.8 

3.8 

0.4 

1.36 

3502 

1.266 

3.6 

3.6 

0.4 

1.29 

E-862 

1.2 

2.7 

2.7 

0.4 

1.0 


(a) 



A 

B 

C 

D 



Figure 6. — FE mesh and boundary conditions for the (a) axial and (b) transverse coupons. 


(particularly for in-plane loading). The process of determining the UD lamina properties was repeated 
depending on the number of UD fiber volume fractions present in the discretization. For example, the 
absorbed matrix model would involve computing effective properties using MAC/GMC for 37.5, 73, and 
80 percent fiber volume fractions. The layered shell formulation within LS-DYNA would use laminate 
theory to determine the overall subcell response based upon the provided UD properties. This bottom-up 
process was used to compute the UD ply properties for four different material systems (T700 Toray fibers 
with different matrix systems: Cytec PR520, Cytec 5208, Flexcel 3502, Epon Epikote Resin 862, hereby 
referred to as E-862). 


3.2 Finite Element Model 

This section highlights the implementation of the layered subcell approach for modeling the braided 
composite within the FE code, LS-DYNA. Each of the individual subcells were modeled as unique parts 
composed of a layered, Belytschko-Tsay quadrilateral shell using the *INTEGRATION_SHELL 
keyword. Each layer was assigned a unique integration point and associated thickness according to the 
parameters provided in Table 2. Tensile coupon tests were simulated in both the axial direction (parallel 
to the axial tows) and the transverse direction (perpendicular to the axial tows). Axial coupon specimen 
geometries were 8 by 1.412 by 0.125 in. (H by W by D, w/240 elements) whereas the transverse coupon 
geometries were 12.002 by 1.4 by 0.125 in. (476 elements). These dimensions corresponded to a 6-ply 
specimen of the braided composite. Each shell element thus represented 6 plies of the braided composite 
and the total thic kn ess of the shell was set to 0. 125 in. It should be noted that in the real composite braid, 
nesting of axial fibers occurs and results in imperfect stacking. These effects are currently neglected and 
perfect layer stacking is assumed. The finite element models for the tensile tests are provided in Figure 6. 
The compression specimen geometries were 1 by 1.412 by 0.125 in. and 1.4 by 0.8714 by 0.125 in. 

(H by W by D) for the axial and transverse compression tests, respectively. The meshes for the axial and 
transverse compression tests used to compute the compressive moduli are shown by highlighted boxes in 
Figure 6. A white arrow in Figure 6, and in figures to follow, identifies the direction of the axial fiber tows. 
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4.0 Coupon Test Results 

The finite element results from the simulated coupon tests were used to compute effective properties 
of the braided composite for various matrix systems. Cross-sectional forces acting at the center gage of 
the specimen are used to compute the uniaxial stresses, whereas crosshead displacements are used for 
strain calculations. The stress-strain curves are then used to compute the moduli for the axial and 
transverse directions in both tension and compression. Table 5 displays the results for the E-862 and 3502 
material systems. Shown are the predicted moduli for the pure matrix model, absorbed matrix model, 
previous subcell models (where available), and experimental test data. Table 6 presents similar data for 
the PR520 and 5208 material systems. General trends in the moduli predictions were that the axial tensile 
stiffnesses were slightly over-estimated (~2.4 percent) and the transverse tensile stiffnesses were 
generally under-predicted (~18 percent). Despite little variation in the experimental transverse moduli 
across the material systems, the simulations were sensitive to the matrix modulus. The pure matrix model 
was slightly stiffer in all cases, due to the localization of the fibers in the 80 percent UD plies and its 
inherent effect on the subcell response. Both investigated subcell models adequately captured the stiffness 
properties even when compared against previous subcell works listed for the E-862 and PR520 systems in 
Tables 5 and 6. The results for the compressive moduli varied across the material systems and between 
the transverse and axial directions. The transverse compression moduli were generally in better agreement 
(~6 percent), however, the axial compression moduli were only slightly over-predictive (~13 percent). For 
the case of compression, a similar trend in matrix moduli sensitivity was found for the transverse case, 
and an apparent softening in the absorbed matrix model in comparison to the pure matrix model was 
noticeable. Considering the purely predictive, bottom-up method of the proposed subcell models the moduli 
predictions were well within a reasonable range. Both the pure matrix and absorbed matrix models are 
investigated for their out-of-plane properties and verified against single ply coupon tests in the next section. 


TABLE 5.— COMPARISONS OF COMPUTED MODULI WITH EXPERIMENTAL DATA FOR 
THE E-862 AND HEXCEL 3502 MATERIAL SYSTEMS 



E-862 


3502 j 

Pure 

matrix 

Absorbed 

matrix 

Previous works, 
(Ref. 30) 

Test, 
(Ref. 34) 

Pure 

matrix 

Absorbed 

matrix 

Test, 
(Ref. 30) 

EllT(GPa) 

48.92 

47.32 

51.10 

46.9 


49.52 

48.41 

47.2 

E22T (GPa) 

32.36 

31.43 

38.60 

41.6 


34.53 

32.39 

40.6 

EllC(GPa) 

50.31 

48.64 


41.4 


51.45 

49.26 

41.8 

E22C (GPa) 

38.44 

34.80 


42.7 


40.94 

36.86 

41.1 


TABLE 6.— COMPARISONS OF COMPUTED MODULI WITH EXPERIMENTAL DATA FOR 
THE CYTEC 5208 AND PR520 MATERIAL SYSTEMS 



5208 


PR520 ! 

Pure matrix 

Absorbed 

matrix 

Test 

(Ref. 34) 

Pure 

matrix 

Absorbed 

matrix 

Previous works 
(Refs. 30 and 27) 

Test 

(Ref. 34) 

EllT(GPa) 

49.858 

48.905 

47.5 


51.496 

49.22 

51.1/44.9 

47.7 

E22T (GPa) 

35.351 

35.377 

41.3 


36.667 

36.81 

38.7/51.2 

42.7 

El 1C (GPa) 

51.744 

50.037 

44.6 


52.154 

51.20 


47.6 

E22C (GPa) 

41.42 

39.790 

38.6 


42.434 

41.76 


42.0 
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5.0 Single Ply Coupon Tests and Simulations 

Previous works have investigated the out-of-plane displacements in single-ply coupon tests of the 
braided composite under axial and transverse tension (Ref. 2) using Digital Image Correlation (DIC) on 
the T700/E862 braided composite. The displacements were a result of the locally unsymmetric nature of 
the braided composite. The DIC images of the out-of-plane displacements are reproduced in Figures 7(a) 
and (b). In light of these results, it was conjectured that a modeling approach best capturing this phenome- 
non would be advantageous for two reasons. One, it would provide a means of quantifying the flexural 
property predictions of the subcell model. Second, it would suggest that the subcell model may be 
sufficient in capturing the effects of inter-ply delamination hypothesized to be caused by these out-of- 
plane deformations. These edge-related delaminations drastically affect the stress/strain distribution in the 
gage section of standard test coupons (Ref. 2) and accurate simulations of these phenomena may provide 
insight into obtaining valuable strength data from material characterization tests. 


5.1 Subcell Flexural Properties 

CLPT was used initially to determine the anticipated response of the two subcell discretizations (pure 
matrix and absorbed matrix models). This analysis would provide insight into the out-of-plane properties 
of the subcells and would allow for meaningful interpretation of the single ply coupon simulations. Two 
subcell characteristics were identified as potential contributors to the out-of-plane deformations and 
studied using CLPT: the tension/twist coupling of subcells A and C and the flexural modulus of subcell B 
and D. To quantify the tension twist coupling under uniaxial loading, relations (1) and (2) below were 
used- whereby assumptions of uniaxial stress and zero membrane moments provided a relation between 
membrane strain, s x or e v , and the twist curvature, k xy . The x-direction corresponds to that of the axial tows 
(0°), and the y-direction corresponds to the transverse direction perpendicular to the axial tows. 

Equation (1) represents the tension-twist coupling for uniaxial tension in the x-direction. Conversely, 
Equation (2) is used for the case of uniaxial tension in the y-direction. The quantities T x and T y represent 
the linear relation between membrane strain and curvature for the cases of uniaxial loading in the axial 
and transverse directions, respectively. The quantities B i6 and B 2 6 are the only non-zero components of 
the extensional-bending coupling matrix of the subcell laminate; it is for this reason the curvatures are 
reduced to only k xv for the case of uniaxial tension. The computed results for the two discretization 
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Figure 7. — z+ displacement contours obtained using DIC for (a) single ply 
transverse tensile tests and (b) single ply axial tensile tests. 


TABLE 7.— TENSION TWIST COUPLING 


IN SUBCELL A 


Model 

r* 

(l/m) 

T* 

(l/m) 

Pure matrix 

-5800.64 

-1216.6 

Absorbed matrix 

-5983.83 

-1208.43 


TABLE 8.— PREDICTED FLEXURAL MODULI 
OF SUBCELL B 


Model 

E J X 

(GPa) 

E i 

r, y 

(GPa) 

E * 

(GPa) 

Pure matrix 

9.48 

21.37 

7.91 

Absorbed matrix 

6.53 

21.01 

10.4 


TABLE 9.— COMPARISON OF BENDING 
STIFFNESS MATRICES FOR SUBCELL B 


Pure matrix model 

Absorbed matrix model 

0.2324 

0.2921 

-0.1558 

0.2966 

0.1403 

0.5667 

0.2921 

0.9500 

-0.4372 

0.1403 

0.2112 

0.4890 

-0.1558 

-0.4372 

0.3201 

0.5667 

0.4890 

0.2203 


models, shown in Table 7, highlight the minute differences between the pure and absorbed matrix models 
in terms of the quantities T x and T y . The flexural moduli of subcell B for the two discretizations were 
computed using Equation (3) and the results are presented in Table 8. The flexural moduli in the x and xy 
directions varied considerably between the subcell discretizations, however the flexural moduli for the 
transverse, y, direction were nearly identical for both models. This similarity in transverse flexural moduli 
was contrary to the unambiguously different D matrices (bending-stiffness matrices) computed for the 
two methodologies, shown in Table 9. 


5.2 Numerical Verification of Single Subcell Responses 

The curvatures seen in Figure 7(a) and (b) were hypothesized to be a direct consequence of the 
tension-twist coupling present in the unsymmetric subcells A and C. Thus, it was necessary to ensure the 
LS-DYNA implementation of the subcell models captured this phenomenon and to determine any mesh 
and/or element-type dependencies. Due to little expected differences between the two subcell 
discretizations in the tension/twist responses of subcells A and C, only the pure matrix model was used in 
the simulations to follow and shown in Figures 8 and 9. Figure 8(a) and (b) present the contour plots of 
out-of-plane displacement for subcell A under uniaxial tension. As can be seen from the contours, the 
displacements in the y and z directions were unrestrained. There was noticeable smoothing in the contours 
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Figure 8. — Contour plots of subcell A’s out-of-plane z-displacement due to 
uniaxial tension in the axial (x) direction for (a) a coarse mesh and (b) a 
more refined mesh. Displacements in the y and z directions are unrestrained 
on all boundaries. The red and blue spectrums on the contour plots represent 
positive and negative z-displacements, respectively. Results are qualitative, 
thus no fringe scales are displayed. Shell elements are Belytschko-Tsay. 



Figure 9. — Contour plots of subcell C’s out-of-plane z-displacement due to uniaxial tension in the axial (x) direction. 
The four contour plots show comparisons for unrestrained y and z displacements with (a) Belytschko-Tsay and 

(b) Belytschko-Leviathan shell elements, and fixed y and z displacements at the top and bottom boundaries for 

(c) Belytschko-Tsay and (d) Belytschko-Leviathan shell elements. 


as a function of increased mesh dependency; however, the size and magnitude of the deformation were 
unchanged. The out-of-plane deformations were as expected from the laminate analysis reported 
previously with only twist curvature, k„, present in the case of uniaxial tension. The various contour plots 
in Figure 9 for subcell C highlight the effects of edge boundary conditions on the displacements as well as 
the influence of shell element formulation. Two shell formulations were investigated, Belytschko-Tsay 
and Belytschko-Leviathan, based on recommendations in literature for modeling unbalanced and 
unsymmetric lamina in LS-DYNA (Ref. 35). In the case of free y and z displacements on the top and 
bottom surfaces of the subcell, the two element formulations produced nearly identical results. When the 
y and z displacements were fixed on the top and bottom boundaries, a reversal in twist was observed. Due 
to the complex nature of the induced curvatures, only Belytschko-Leviathan shells obtained smooth 
contour plots seen in Figure 9(d). In contrast, the Belytschko-Tsay elements had noticeable sharp, lined 
features which appear to be a function of the element size, highlighted in Figure 9(c). It should also be 
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mentioned that comparing Figure 9(a) with Figure 8(a) indicates that the out-of-plane deformation of 
subcells A and C have reflective symmetry with respect to each other. Based on the performance of the 
Belytschko-Leviathan shells, it was deemed that the LS-DYNA layered shell implementation of the 
subcell discretizations were sufficient to capture the necessary tension-twist mechanisms. Belytschko- 
Leviathan shell elements are utilized in the single ply coupon simulations to follow. 


5.3 Single Ply Coupon Simulations 

Full coupon simulations were performed for both axial and transverse tension using similar coupon 
dimensions from the earlier tensile specimens. In this study, the thicknesses of the shells were reduced to 
that of a single ply (-0.022 in.). The coupon mesh was refined to include 4 elements per subcell, shown in 
Figure 10(b), in comparison to that of the standard coupon tests in Figure 6(a) and (b). Belytschko 
Leviathan shell elements are used to better account for the tension-twist coupling of the unsymmetric 
subcells as discussed in the previous section. Both coupon specimens are loaded to 1 percent strain for 
which the composite response in the experiments was purely elastic. Thus, the currently elastic material 
models for the UD plies in each subcell are again deemed appropriate for capturing the out-of-plane effects. 

Transverse coupon results are displayed in Figure 10(a) and (c) for the pure matrix model and can be 
compared to mesh/subcell location using Figure 10(b). The two contour plots in Figure 10(a) and (c) 
represent two alternating stable positions in the predicted solution, which tend to alternate through 
loading. Figure 10(a) contained the major curvatures and displayed similarities to the DIC results in 
Figure 7(a). Minor curvatures of opposite sign seen in Figure 10(c) also existed between the major 
curvatures that could be observed in the DIC data. Two distinct remarks could be made, qualitatively, 


Axial 

direction 

L_ Loading direction 


; 


A 

B 

C 

D 





Figure 10. — Contour plots of z-displacement for a transverse coupon specimen loaded in the y-direction. Axial fiber 
tows are oriented in the x-direction. (a) Stable position with only major curvatures displayed, (b) coupon mesh for 
comparing location of displacements with respect to the subcells, and (c) a stable position with both minor and 
major curvatures. The meshes in (a) through (c) are aligned for comparisons. Results are qualitative, thus no 
fringe scales are displayed. 
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Figure 1 1 . — Contour plots of z-displacement for a transverse coupon specimen loaded in the y-direction and axial 
fiber tows oriented in the x-direction. The B/D subcells closest to left and right boundaries are labeled accordingly. 
The dashed box highlights a transitional region where the location of major curvatures is switching. 

from the contour results. First, edge twisting resulted from a collective effect of subcells A and C and the 
reflective symmetry of their tension-twist coupling. Consequently, significant edge displacements, 
although rooting from the twisting of subcells A and C, were manifested in subcells B and D. In 
Figure 10(c), the major curvatures occurred in the regions of subcell B, whereas the minor curvatures 
occurred within subcells D. This feature was determined to be a consequence of boundary conditions and 
the relative adjacency of either of these subcells with the restrained edge. For example, modifying the 
order of subcells so that subcell C would be the first from the constrained edge — that is, if looking left to 
right in Figure 10(b) — would move the location of the major curvatures to subcells D. This shifting of 
subcells was performed and the single ply transverse tension simulation is given in Figure 11. The 
z-displacement contours shown in Figure 1 1 provide an example of competing major/minor curvatures 
between subcell B and D by varying the subcell closest to the restrained edge (labeled on the figure). The 
left hand side of Figure 1 1 features dominant positive (red) displacement curvatures at the top of the 
model, whereas beyond the dashed transitional region marked on the figure the dominant displacement 
curvatures at the top of the model are negative (blue). Though the source of major/minor curvatures was 
found to be a function of boundary condition, the alternating stable contours between those in 
Figure 10(a) and (c) were thought to be a construct of the numerical implementation. The contours seen 
experimentally were entirely stable and matched closely to that of Figure 10(c), and thus supported this 
deduction. These results highlight the importance of the finite element model’s constraints on the out-of- 
plane displacement contours in the modeling approach, a phenomenon which must be acknowledged 
when comparing experimental and numerical results. 

A quantitative comparison of the z-displacements for the two subcell models is displayed in 
Figures 12(a) to (c) with the full field DIC data presented again in Figure 12(d). The results are presented 
using the stable contour shape from Figure 10(c), which are found to match closely to those seen 
experimentally. Two different element sizes are provided for the absorbed matrix model, Figure 12(a) and 
(b), to show the effect of a mesh refinement on smoothing the displacement contours. Note that all four 
contour plots match to the fringe levels shown at the right of Figure 12. As expected from the single 
subcell tests, the finer discretization, although smoothing the contours, did not alter the magnitude of the 
z-displacements. Comparing the two subcell models, the pure matrix model showed a significant 
reduction in out-of-plane deformation. Despite the pure matrix model showing qualitatively accurate 
curvatures in Figures 10 and 11, the contours were undeipredicted in comparison to the experimental 
results. Based on the earlier investigations using CLT, the probable cause of the suppression of the out-of- 
plane displacements in the pure matrix model was the divergent properties of the bending stiffness matrix 
components between the two discretizations. The pure matrix model had a 145 percent higher flexural 
modulus in the axial fiber direction than the absorbed matrix counterpart, despite there being almost no 
differences in the tension-twist coupling of subcell A between the two models. Thus, the stiff outer UD plies 
in subcell B of the pure matrix model are acting like composite sandwich skin layers, suppressing the out- 
of-plane deformation. The z-displacements predicted by the absorbed matrix model for both element 
discretizations in Figure 12(a) and (b) were in good agreement with the experimental results in Figure 12(d). 
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Figure 12. — Comparisons of z-displacement for absorbed matrix model with both 
(a) a coarse mesh and (b) a fine mesh; (c) pure matrix model with a coarse mesh; 
and (d) experimental DIC data (Ref. 2). The units of the LS-DYNA contour legend 
are in inches and are equivalent to that of the DIC results (originally presented in 
Figure 7(a) in mm). 


The z-displacement contour plots for axial coupon simulations are shown in Figure 13. The axial 
coupon was meshed coarsely with two elements per subcell width, similar to that of the transverse case 
shown in Figure 12(a). Although there were out-of-plane displacements occurring near the gripped ends 
of the coupon, readily apparent in the absorbed matrix model — the results did not mimic the convex/ 
concave undulations throughout the body of the coupon as seen experimentally in Figure 7(b). The nature 
of the tension-twist coupling of subcells A and C permitted twist deformations only at the comers of the 
subcells due to anti- symmetry, thus almost no out-of-plane deformations were observed in the interior of 
the model. It must be mentioned that the contour plots from the LS-DYNA simulations in Figure 13(a) 
and (b) are on a scale l/6th that of the contour plots in the experimental data in Figure 7(b). Thus, 
although these deformations are present within the single ply and may contribute to internal stresses 
within the 6-layer braid, the effects of the axial undulations are not as significant to coupon damage, 
unlike the edge twist in the transverse single ply coupons. Consequently, it was concluded that the 
mechanisms driving the concave/convex curvatures seen experimentally in the axial coupon tests arose 
from physical mechanisms other than the locally unsymmetric nature of the single-ply braided composite. 
It is thought that capturing concave/convex curvatures in the axial coupon simulations requires a 
methodology preserving fiber continuity. Due to the nature of the discretization and the resulting 
discontinuity of fiber tows, capturing curvature phenomena like those seen in Figure 7(b) is a limitation to 
the current subcell approach. 
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Figure 13. — Contour plots of z-displacement for axial coupon tests 
(load applied in fiber, or x-direction) for (a) the absorbed matrix 
model and (b) the pure matrix model. Units of the contour legend 
are in inches. 


6.0 Conclusion 

Based on the results from the single ply coupon tests, the absorbed matrix model was found to be an 
excellent candidate for future advancements in the semi-analytical approach. The adequate moduli 
prediction in both tension and compression, as well excellent agreement in out-of-plane tension-twist 
displacements, under global transverse loading, were based purely on constituent properties, thus 
highlighting the efficacy of the model to be extended to other material systems. The analysis of the 
flexural properties of the subcell discretizations was a novel contribution and provided insight into the 
mechanisms driving the out-of-plane deformations seen in the transverse coupon tests. It was found that 
the unsymmetric A and C subcell laminates were the direct cause of the unique material response. The 
excellent correlation with the experimental DIC data for the transverse coupon case supported the subcell 
discretization utilized in the absorbed matrix model, which left subcells B and D as symmetric laminates. 
The single ply coupon results also verified the assumption that the matrix was best modeled in series with 
surrounding braider tows, rather than in parallel as in the pure matrix model. The subcell model, however, 
was not able to capture the axial out-of-plane responses in the single ply simulations and was attributed to 
a limitation of the overall discretization scheme. Although very slight differences were found between the 
two subcell models for in-plane properties, significant differences in the axial flexural modulus 
contributed to varying tension-twist responses in the single ply. Lastly, the current method which preserves 
the through-thickness subcell properties should be beneficial when performing impact simulations. 


NASA/TM— 2013-217875 


15 




The current semi-analytical subcell approach has several advantages over that of the previous 
iterations, summarized in Section 2.0. Most important, however, is the possibility of using the current 
subcell implementation to better understand the braided coupon’s edge related damage, currently a 
limiting constraint in determining material strength parameters. In conjunction with understanding coupon 
damage mechanisms, future work will look to develop a top-down methodology in order to determine the 
necessary UD ply strengths for the absorbed matrix subcell discretization. 
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Appendix — Straight Line Braided Composite Model 

The straight line model (Ref. 29) for determining the necessary UD ply thicknesses and computing 
the volume of fibers in the axial and braider directions contained in the subcells is presented here. The 
subcell widths W A and Wg are measured as shown in Figure Al. Figure A1 also shows the amounts of 
+60° braider tow ( dashed yellow box) approximated to exist inside of the subcell dimensions highlighted 
by the green box. 

Equations (4) to (6) are used to compute the volume of fiber contained in subcells A or B (the 
quantity of fibers in C and D will be identically equivalent). The subscripts a and b designate axial or 
braider tow and A and B subscripts designate the corresponding subcell. Superscript / represents volume 
of fibers whereas the subscript / represents volume fraction. The tow subscript represents the fiber tow, m 
the matrix layer, L the length of the subcell shown in Figure Al, and h the height of the subcell through 
the thickness as shown in Figure A2. 
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Equations (7) and (8) compute the UD ply thicknesses for the pure matrix model. In both subcell 
models, the volume fractions of the three different UD plies (axial tow in A, braider tows in A, and 
braider tows in B) are computed according to Equations (9) to (1 1). The thickness of the axial UD in 
subcell A is assumed to be equivalent to the original pure matrix model; however, the braider UD plies in 
both A and B are modified according to Equations (12) and (13). 
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Figure A1. — Two unit cells of the braided composite are shown with the widths of subcells A and B 
labeled. The green boxes represent the size of a subcell, whilst the dotted yellow boxes represent 
the amount of braider tow contained within the subcell volume as a straight section of fiber tow. 



Figure A2. — Schematic of the braider fiber tow lengths approximated by the straight line mode 
for subcells A and B. 
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